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Abstract. This paper presents a technique for reducing speckle in Po- 
larimetric Synthetic Aperture Radar (PolSAR) imagery using Nonlocal 
Means and a statistical test based on stochastic divergences. The main 
I I objective is to select homogeneous pixels in the filtering area through 

f-H statistical tests between distributions. This proposal uses the complex 

^^ Wishart model to describe PolSAR data, but the technique can be ex- 

C^ tended to other models. The weights of the location- variant linear filter 

, ^, are function of the p- values of tests which verify the hypothesis that two 

samples come from the same distribution and, therefore, can be used 
^H to compute a local mean. The test stems from the family of (/i-(/)) di- 

^ vergences which originated in Information Theory. This novel technique 

■^T was compared with the Boxcar, Refined Lee and ID AN filters. Image 

C^ quality assessment methods on simulated and real data are employed to 

^P validate the performance of this approach. We show that the proposed 

filter also enhances the polarimetric entropy and preserves the scattering 



■^^ information of the targets. 

o 

CO Keywords: Hypothesis testing, Information theory, Multiplicative noise, 

^"^ PolSAR imagery. Speckle reduction. Stochastic distances. Synthetic Aper- 

^ ture Radar 



^ 1 Introduction 



Among the remote sensing technologies, Polarimetric Synthetic Aperture Radar 
(PolSAR) has achieved a prominent position. PolSAR imaging is a well-developed 
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coherent microwave remote sensing technique for providing large-scale two-di- 
mensional (2-D) high spatial resolution images of the Earth's surface dielectric 
properties [20 . 

In SAR systems, the value at each pixel is a complex number: the amplitude 
and phase information of the returned signal. Full PolSAR data is comprised of 
four complex channels which result from the combination of the horizontal and 
vertical transmission modes, and horizontal and vertical reception modes. 

The speckle phenomenon in SAR data hinders the interpretation these data 
and reduces the accuracy of segmentation, classification and analysis of objects 
contained within the image. Therefore, reducing the noise effect is an important 
task, and multilook processing is often used for this purpose in single- and full- 
channel data. In the latter, such processing yields a covariance matrix in each 
pixel, but further noise reduction is frequently needed. 

According to Lee and Pottier [20 , Polarimetric SAR image smoothing re- 
quires preserving the target polarimetric signature. Such requirement can be 
posed as: (i) each element of the image should be filtered in a similar way to 
multilook processing by averaging the covariance matrix of neighboring pixels; 
and (ii) homogeneous regions in the neighborhood should be adaptively selected 
to preserve resolution, edges and the image quality. The second requirement, 
i.e. selecting homogeneous areas given similarity criterion, is a common problem 
in pattern recognition. It boils down to identifying observations from different 
stationary stochastic processes. 

Usually, the Boxcar filter is the standard choice because of its simple design. 
However, it has poor performance since it does not discriminate different targets. 
Lee et al. |17|18 propose techniques for speckle reduction based on the multi- 
plicative noise model using the minimum mean-square error (MMSE) criterion. 
Lee et al. [19] proposed a methodology for selecting neighboring pixels with sim- 
ilar scattering characteristics, known as Refined Lee filter. Other techniques use 
the local linear minimum mean-squared error (LLMMSE) criterion proposed by 
Vasile et al. [37 , in a similar adaptive technique, but the decision to select homo- 
geneous areas is based on the intensity information of the polarimetric coherency 
matrices, namely intensity-driven adaptive- neighborhood (ID AN). 

Qetin and Karl [4 presented a technique for image formation based on regu- 
larized image reconstruction. This approach employs a tomographic model which 
allows the incorporation of prior information about, among other features, the 
sensor. The resulting images have many desirable properties, reduced speckled 
among them. Our approach deals with data already produced and, thus, does 
not require interfering in the processing protocol of the data. 

Osher et al. [25] presented a novel iterative regularization method for inverse 
problems based on the use of Bregman distances using a total variation denoising 
technique tailored to additive noise. The authors also propose a generalization for 
multiplicative noise, but no results with this kind of contamination are show. The 
main contributions were the rigorous convergence results and effective stopping 
criteria for the general procedure, that provides information on how to obtain 
an approximation of the noise-free image intensity. Goldstein and Osher [15] 
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presented an improvement of this work using the class of Li -regularized opti- 
mization problems, that originated in functional analysis for finding extrema of 
convex functionals. The authors apply this technique to the Rudin-Osher-Fatemi 
model for image denoising and to a compressed sensing problem that arises in 
magnetic resonance imaging. Our work deals with full polarimetric data, for 
which, to the best of our knowledge, there are no similar results that take into 
account its particular nature: the pixels values are definite positive Hermit ian 
complex matrices. 

Soccorsi et al. [28 presented a despeckling technique for single-look com- 
plex SAR image using nonquadratic regularization. They use an image model, a 
gradient, and a prior model, to compute the objective function. We employ the 
full polarimetric information provided by the multilook scaled complex Wishart 
distribution. 

Chambolle [5 proposed a Total Variation approach for a number of problems 
in image restoration (denoising, zooming and mean curvature motion), but under 
the Gaussian additive noise assumption. 

Li et al. [21 propose the use of a particle swarm optimization algorithm and 
an extension of the curvelet transform for speckle reduction. They employ the 
homomorophic transformation, so their technique can be used either in ampli- 
tude or intensity data, but not in complex-valued imagery, as is the case we 
present here. 

Wong and Fieguth [41] presented a novel approach for performing blind 
decorrelation of SAR data. They use a similarity technique between patches 
of the point-spread function using a Bayesian least squares estimation approach 
based on a Fisher-Tippett log-scatter model. In a similar way, S0lbo and Eltoft [29] 
assume a Gamma distribution in a wavelet-based speckle reduction procedure, 
and they estimate all the parameters locally without imposing a fixed number 
of looks (which they call "degree of heterogeneity") for the whole image. 

Buades et al. [3j proposed a methodology, termed Nonlocal Means (NL- 
means), which consists in using similarities between patches as the weights 
of a mean filter; it is known to be well suited for combating additive Gaus- 
sian noise. Deledalle et al. [TT] applied this methodology to PolSAR data using 
the Kullback-Leibler distance between two zero-mean complex circular Gaus- 
sian laws. Following the same strategy, Chen et al. [6] used the test for equality 
between two complex Wishart matrices proposed by Conradsen et al. ^. 

This paper proposes a new approach for speckle noise filtering in PolSAR 
imagery: an adaptive nonlinear extension of the NL-means algorithm. This is 
an extension of previous works |33|34j . where we used an approach similar to 
that of Nagao and Matsuyma [23]. Overlapping samples are compared based 
on stochastic distances between distributions, and the p-values resulting from 
such comparisons are used to build the weights of an adaptive linear filter. 
The goodness-of-fit tests are derived from the divergences discussed by Frery 
et al. [M] and Nascimento et al. [24]. The new proposal is called Stochastic 
Distances Nonlocal Means (SDNLM), and amounts to using those observations 
which are not rejected by a test seeking for a strong stationary process. 
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This paper is organized as follows. First, we summarize the basic principles 
that lead to the complex Wishart model for full polarimetric data. In Section|3]we 
recall the Nonlocal Means method. Our approach for reducing speckle in PolSAR 
data using stochastic distances between two complex Wishart distributions is 
proposed in Section[4J Image Quality Assessment is briefly discussed in SectionjS] 
Results are presented in Section |6) while Section [7| concludes the paper. 

2 The Complex Wishart Distribution 

PolSAR imaging results in a complex scattering matrix, which includes intensity 
and relative phase data [14]. Such matrices have usually four distinct complex 
elements, namely Syy^ Svhi Shvi ^ind Shhi where H and V refer to the horizon- 
tal and vertical wave polarization states, respectively. In a reciprocal medium, 
which is most common situation in remote sensing, Svh = Shv so the com- 
plex signal backscattered from each resolution cell can be characterized by a 
scattering vector Y with three complex elements [36 . 

Thus, we have a scattering complex random vector Y = [Shh^Svh^ SyvY^ 
where [■]* indicates vector transposition. In general, PolSAR data are locally 
modeled by a multivariate zero-mean complex circular Gaussian distribution 
that characterize the scene reflectivity [35 36 , whose probability density function 
is 

where | • | is the determinant, and the superscript '*' denotes the complex con- 
jugate of a vector; U is the covariance matrix of Y. This distribution is defined 
on C^. The covariance matrix U, besides being Hermitian and positive definite, 
has all the information which characterizes the scene under analysis. 

Multilook processing enhances the signal-to-noise ratio. It is performed av- 
eraging over L ideally independent looks of the same scene, and it yields the 
sample covariance matrix Z given, in each pixel, by Z = L~^ X]/^=i '^i^* •> where 
L is the number of looks. 

Goodman [16] proved that Z follows a scaled multilook complex Wishart 
distribution, denoted by Z ^ W(^7 L)^ and characterized by the following prob- 
ability density function: 

r3L|^/|L-3 

/z(Z';i;,L) = ^l-^exp{-Ltr(i;-iZ')}, (1) 

where, for L > 3, r^{L) = tt^ Y[i=o^i^ ~ 0^ ^(') ^^ ^^^ gamma function, tr(-) 
is the trace operator, and the covariance matrix Z is given by 



^ = E{rr**} = 



EIShhS^jj} EIShhSyjj} E{SHHSln/} 
E{SvhS^jj} E{SvhSyjj} E{SvhS^} 
E{SwS*jji^} E{SwSyH} E{SvvSlrY} 



where E{-} denote expectation. Anfinsen et al. [2 removed the restriction L > 3. 
The resulting distribution has the same form as in ([l]) and is termed the "relaxed" 
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Wishart. We assume this last model, and we allow variations of L along the 
image. 

The support of this distribution is the cone of positive definite Hermitian 
complex matrices p!2l. 

The parameters are usually estimated by maximum likelihood (ML) due to 
its statistical properties. Let Zr be a random matrix which follows a W(X',L) 
law. Its log-likelihood function is given by 

£r{^,L) = 3L\ogL^{L-3)\og\Zr\ -Llogl^"! -31og7r- 

2 

Y,^ogr{L-q)-LtT{U-'Zr), 

q=0 

resulting in the following score function: 

LveciU-^ZrlJ-^ - X-i) 



^^^ ' 3(logL + l)+log|Z,|-log|^|-E-=oV^'(^-j')-tr(^-'^r);' 

where i/j^ is the digamma function. Let {Zi, Z2, . . . , Ztv} be an i.i.d. random 
sample of size N from the W(^, L) law. The ML estimator (X", L) of its param- 
eters is IJ = Z = N~^U^^iZr and the solution of 

N 2 

31og£+-^log|Z,|-log|Z|-^^0(i-9)=0. (2) 

The case L < q was also treated by Anfinsen et al. [2]. 

3 Nonlocal Means 

The NL-means method was proposed by Buades et al. [3 based on the redun- 
dancy of neighboring patches in images. In this method, the noise-free estimated 
value of a pixel is defined as a weighted mean of pixels in a certain region. Un- 
der the additive white Gaussian noise assumption, these weights are calculated 
based on Euclidean distances which are used to measure the similarity between 
a central region patch and neighboring patches in a search window. The filtered 
pixel is computed as: 

9{x^y) = F^ 7 — ^ , (3) 

where w{u^v) are the weights defined on the search window W. The resulting 
image g is, thus, the convolution of the input image / with the mask w — 
w^ / Y^ w'iu^ v). The factor w'{u^ v) is inversely proportional to a distance between 
the patches, and is given by 

w'{u,v) = exp{-i ^|/(p„(fc)) - f{p,{k))\'}, 
keP 
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where h > controls the intensity, in a similar way the temperature controls 
the Simulated Annealing algorithm, f{pu{k)) and f{pv{k)) are the observations 
in the k-th pixels of the patches centered in u and v, respectively. When h -^ 
oo the weights tend to be equal, while when h ^ they tend to zero unless 
f{pu{k)) = f{pv{k)). In the first case the filter becomes a mean over the search 
window; in the last case, the filtered value will remain unaltered. 

The NL-Means proposal and its extensions rely on computing the weights 
of a convolution mask as functions of similarity measures: the closer (in some 
sense) two patches are, the heavier the contribution of the central pixel to the 
filter. 

Deledalle et al. [10 analyze several similarity criteria for data which depart 
from the Gaussian assumption. In particular, the authors consider the Gamma 
and Poisson noises, because they are good image models. Deledalle et al. [9\ 
extended the NL-means method to speckled imagery using statistical inference 
in an iterative procedure. The authors derive the weights using the likelihood 
function of Gaussian and square root of Gamma (termed "Nakagami-Rayleigh" ) 
noises. This idea is extended by Deledalle et al. [11 to PolSAR data under the 
complex Gaussian distribution, and by Su et al. [30 to multitemporal PolSAR 
data. 

Chen et al. [6 presented a NL-Means filter for PolSAR data under the com- 
plex Wishart distribution. The authors employ the likelihood ratio test of equal- 
ity between two W(X', L) laws with the same number of looks L, as presented 
by Conradsen et al. [8 , to calculate the weights between the patches. 

Our proposal addresses the problem in a more general fashion: we use the 
p- value of a goodness-of-fit test between two samples. The tests are derived using 
an Information Theory approach to compute stochastic divergences, which are 
turned into distances and then scaled to exhibit good asymptotic properties: 
they obey a x^ distribution. The weights are computed with a soft threshold 
which incorporates all the observations which were not rejected by the test, and 
some of the others. As presented in the next section, all these elements can be 
easily generalized to other models, tests and weight functions. 

4 Stochastic Distances Filter 

4.1 The weights 

In this paper we use a 7 x 7 search window. The shape and size of the neighboring 
patches and the central patch are the same: squares of 3 x 3 pixels. The central 
patch, with center pixel Zi, is thus compared with 24 neighboring patches, whose 
center pixels are Z^, z = 2, . . . , 25, as illustrated in Figure [l] 

The estimate of the noise-free observation at Zi is a weighted sum of the 
observations at Z2, . . . , Z25, being each weight a function of the p- value (p(l, i)) 
observed in a test of same distribution between two complex Wishart laws: 

if p{l,i) > T], 
^p(l,i)-lif f <p(l,z)<7?, (4) 

otherwise, 
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Fig. 1. Central pixel Z\ and its neighboring Z'^, % — {2, . . . , 25} with patches 3x3 
pixels and filtering window 5x5 pixels. 



where 77 is the significance of the test, specified by the user. By definition, 
p(l, 1) = 1. This function is illustrated in Figure |2j 



'^(l,0>s 




7y/2 ^1 ^(I'O 

Fig. 2. Weight function. 



In this way we employ a soft threshold instead of an accept-reject decision. 
This allows the use of more evidence than with a binary decision, as the one 
used in |32|34j which was 1 if the sample was not rejected and otherwise. 

When all weights are computed, the mask of convolution coefficients is scaled 
to add one. 

This setup can be generalized in many ways, among them: 



The support, i.e., the set of positions i do not need to be local; it can extend 

arbitrarily. 

The shape and size of the local windows. 
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— The weight function; we opted for the piecewise hnear function presented in 
equation Q as a good compromise between generahty and computational 
cost. 

— The test which produces each p- value. 

The rationale for using the weight function specified in equation Q is the 
following. If the p- values were used instead, the presence of a single sample with 
excellent match to the central sample would dominate the weights in the mask, 
forcing other samples that were not rejected by the test to be practically dis- 
carded. For instance, consider the case p(l,ii) = 0.89, p(l,Z2) = 0-05 and all 
other p-values close to zero. Without the weight function, the nonzero weights 
would be, approximately, 0.46 and 0.03, so the second observation would have 
a negligible infiuence on the filtered value while the first, along with the cen- 
tral value, dominate the result. Using the aforementioned expression, the three 
weights would equal 1/3. This increases the smoothing effect without loosing the 
discriminatory ability. 



4.2 The statistical test 

As defined in the previous section, filtering each pixel requires computing a 
number of goodness-of-fit tests between the patch around the central pixel Z\ 
and the patches surrounding pixels Z^, 2 < z < 25. This will be performed using 
tests derived from stochastic distances between samples. 

Denote Q\ the estimated parameter in the central region Zi, and (^2, • • • , ^25) 
the estimated parameters in the remaining areas. To account for possible de- 
partures from the homogeneous model, we estimate — (Yl^L) by maximum 
likelihood. 

The proposal is based on the use of stochastic distances between the patches. 
Consider that Z\ and Z^, 2 < i < 25, are random matrices defined on the 
same probability space, whose distributions are characterized by the densities 
fz-i^{Z';Oi) and fz^iZ'^Oi)^ respectively, where Oi and Oi are parameters. As- 
suming that the densities have the same support given by the cone of Hermitian 
positive definite matrices A, the h-cj) divergence between fz^ and fzi is given 

by 

where h: (0, 00) -^ [0, 00) is a strictly increasing function with h{0) = and 
h\x) > for every x G R, and (j): (0, 00) -^ [0, 00) is a convex function [26] . 
Choices of functions h and (j) result in several divergences. 

Divergences sometimes are not symmetric. A simple solution, described in 
Frery et al. [14 24 13j , is to define a new measure d^ given by 

du,[Zi,Zi) = . 
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Distances, in turn, can be conveniently scaled to present good statistical prop- 
erties that make them suitable as test statistics 26 : 



Ss{0i,6i) 



2mn 



{m^n)h' {{))(!)" {!) 



dAOi.Oi), 



(5) 



where Oi and Oi are maximum likelihood estimators based on samples size m 
and n, respectively. When Oi = 6i, under mild conditions S^(0i^0i) is asymp- 
totically xif distributed, being M the dimension of 6i. Observing <S?(^i,^i), 

the (asymptotic) p-value of the test is p = Pr(xif > S^iOi^Oi))^ and the null 
hypothesis Oi = Oi can be rejected at significance level r] if p < rj] details can 
be seen in the work by Salicru et al. [26 . Since we are using the same central 
sample for 24 nontrivial tests (notice that S^(0i^0i) = 0). 

Frery et al. |14|13j obtained several distances between W(^, L) distributions. 
The test used in this paper was derived from the Hellinger distance, yielding: 



SniOi^Oi) = 



8mn 
m -\- n 

2 

n 

q=0 



2 . 



(Li+L,)/2 



n 



2 



-q) 



^3Li fSLi 



(6) 



'r{L, - q) r{Li - q) 
provided L > 3. If L < 3, the following expression can be used 



Snieue,) 



8mn 



(^ 



i:i| \Si 



it was also derived by Frery et al. [ 14|13j for the case when Li = L2 = L; in 
this last case, the equivalent number of looks L should be estimated using all 
the available data, for instance as presented by Anfinsen et al. [2j. The setup 
presented in Section |4.1| leads to m = n = 9 since squared windows of side 3 are 
used. 

The tests based on the Kullback-Leibler, Bhattacharyya and Renyi of order 
P G {.1, .5, .9} were also used. They produced almost exactly the same results 
as of the Hellinger distance, at the expense of more computational load. 

Although the distribution of the test statistics given in equation ^ is only 
known in the limit when m, n ^ 00 such that m/{m-\-n) ^ A G (0, 1), it has been 
observed that the difference between the asymptotic and empirical distributions 
is negligible in samples as small as the ones considered in this setup [13pl4j . 

The filter obtained using the p-values produced by the test statistic given 
in equation ^ can be applied iteratively. The complex Wishart distribution is 
preserved by convolutions, and since the number of looks is estimated in every 
pairwise comparison, the evolution of the filtered data is always controlled. This 
also holds for any filter derived from h-(p divergences and the setup presented in 
Section O 

Computational information about our implementation is provided in [A) 
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5 Image Quality Assessment 

According to Wang et al. |39j image quality assessment in general, and filter 
performance evaluation in particular, are hard tasks and crucial for most image 
processing applications. 

We assess the filters using simulated and real input data by means of visual 
inspection and quantitative measures. Visual inspection includes edges, small 
features and color preservation. The measures are computed on the HH^ HV 
and W intensity channels, and they amount to (a) the Equivalent Number of 
Looks (ENL) in homogeneous areas; and (b) the Structural Similarity Index 
Measure (SSIM), proposed by Wang et al. [40] for measuring the similarity be- 
tween two images, (c) the Blind/Referenceless Image Spatial QUality Evaluator 
(BRISQUE), proposed by Mittal et al. [22 to a holistic measure of quality on 
no-reference images. 

The equivalent number of looks N can be estimated by maximum likelihood 
solving equation (pi). The ENL although it is a measure of the signal to noise 
ratio, and it must be considered carefully since it will repute highly blurred 
data as excellent: for this measure, the best result is a completely flat (constant) 
image, i.e., an image with no information whatsoever. The Boxcar filter usually 
has a higher ENL than other proposals because it will always use all the pixels in 
the filtering window and the result is a large loss of spatial resolution. Therefore, 
the equivalent number of looks should be used with care and only to assess noise 
reduction and not image quality in general. 

The SSIM index is a structural information measure that represents the 
structure of objects in the scene, regardless the average luminance and con- 
trast. The SSIM index takes into account three factors: (I) correlation between 
edges; (II) brightness distortion; and (III) distortion contrast. Let / and g be the 
original data and the filtered version, respectively, the SSIM index is expressed 
by 

Cov(/,^) + Ci 2/^ + C2 2afag^Cs 



SSIM(/,^) 



aj dg + Ci /' + ^2 ^ C2 ^/ + ^' + ^3 



where / and g are sample means, a? and a^ are the sample variances, Cov(/, ^f) 
is the sample covariance between / and ^, and the constants Ci,C2,C3 are 
intended to stabilize the index. Following [40], we used the following values: 
Ci = (i^ivC)^, C2 = (i^2>C)^ and Cs = ^2/2, where C is the observed dynamic 
range and both i^i, K2 <C 1; we used Ki = 0.01 and K2 = 0.03. 

The SSIM is defined for scalar- valued images and it ranges in the [—1,1] 
interval, being the bigger value observed is the better result. It was computed 
on each intensity channel as the means over squared windows of side 8, and the 
reported value is the mean over the three channels. This index is a particular 
case of the Universal Quality Index proposed by Wang and Bovik [38^ . We also 
applied this last measure to all the images here discussed, and the results were 
in full agreement with those reported by the SSIM index. 

The BRISQUE is a model that operates in the spatial domain and requires 
no-reference image. This image quality evaluator does not compute specific dis- 



Speckle Reduction in Polarimetric SAR Imagery with SD and NL-means 11 

tortions such as ringing, blurring, blocking, or aliasing, but quantifies possible 
losses of "naturalness" in the image. This approach is based on the principle 
that natural images possess certain regular statistical properties that are mea- 
surably modified by the presence of distortions. No transformation to another 
coordinate frame (DFT, DCT, wavelets, etc) is required, distinguishing it from 
previous blind/no-reference approaches. The BRISQUE is defined for scalar- 
valued images and it ranges in the [0, 100] interval, and smaller values indicate 
better results. 

6 Results 

The proposed filter, termed "SDNLM {Stochastic Distances Nonlocal Means) 
filter" was compared with the Refined Lee (Scattering- model-based), ID AN 
(intensity-driven adaptive-neighborhood) and Boxcar filters. These filters act on 
a well defined neighborhood (as our implementation of the SDNLM): a squared 
search window of side of 5 x 5 pixels, and the former is adaptive (as is the 
SDNLM). 

6.1 Simulated data 

Sampling from the Wishart distribution The simulated data was obtained 
mimicking real data from six classes on a phantom, a segmented image of size 
496 X 496. The original data were produced by a polarimetric sensor aboard the 
R99-B Brazilian Air Force aircraft in October 2005 over Campinas, Sao Paulo 
State, Brazil. The sensor operates in the L-band, and produces imagery with 
four nominal looks. The data were simulated in single look, i.e., in the lowest 
signal-to- noise possible configuration. The covariance matrices are reported m\B\ 
Figures [3 (a) | and [3(b)| show, respectively, the phantom and simulated single- 
look PolSAR image with with false color using \Shh\'^ in the Red channel, \Shv\'^ 
in the Green channel, and |5'wP in the Blue channel. The filtered versions with 
the Boxcar, Refined Lee, IDAN and SDNLM filters are shown in Figures 3(c)[ 



|3(d)[ |3(e) and |3(f)[ respectively. The latter was obtained at the 77 = 90% confi- 



dence level. 

The noise reduction is noticeable in all the four filtered images. The speckled 
effect was mitigated without loosing the color balance, but the Boxcar filter 
produces a blurred image. The Refined Lee filter also introduces some blurring, 
but less intense that the Boxcar filter. The IDAN filter also reduces noise, but 
less effectively that the Refined Lee. The SDNLM filter is the one which best 
preserves edges. This is particularly noticeable in the fine strip which appears 
light, almost white, to the upper left region of the image. After applying the 
Boxcar filter it appears wider than it is. The star-shaped bright spot in the 
middle of the blue field looses detail after Boxcar, Refined Lee and IDAN filters 
are applied, while the SDNLM noise reduction technique preserves its details 
well. The Refined Lee filter introduces a noticeable pixellation effect which is 
evident mainly in the brown areas. 
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The Canny detector was applied to the HH band of ah the available data 
presented in Figure [3) and the results are presented in Figure |4] Even when 
applied to the phantom, the Canny detector fails to produce continuous edges in 
a few situations. In particular, it completely misses the small light features and a 
few edges between regions; see Figure [4 ( a) [ Edge detection with this method is, 
in practice, impossible using the original single-look data, cf. Figure [4(b)] wit h the 
other results, although a few large linear features are visible in the clutter of noisy 
edges. As expected, the edges detected in the Boxcar filtered image are smooth, 
but they neither grant continuity nor identify fine details; this information is 
lost by the filter, see Figure 4(c)| The edges detected in the image processed by 



the IDAN filter are only marginally better than those observed in the original, 
unfiltered, data; see Figure [4(e)| Figures [4(d) [ and [4(f) [ are the edges detected on 
the data filtered by the Refined Lee and SDNLM filters, respectively. Although 
they look alike, it is noticeable that the latter preserves better the small details; 
see, for instance, the star-shaped object to the center-right of the image. It 
appears round in the former, while in the latter it is possible to identify minute 
variations. 

The phantom is available, so it is possible to make a quantitative assessment 
of the results. Table [T] presents the assessment of the filters in the three intensity 
channels HH^ HV and W of the data presented in Figure [3] using the Equivalent 
Number of Looks (ENL) and the SSIM index. 



Table 1. Image quality indexes in the images shown in Figure [S] 



Filter 




ENL 


SSIM Index 




HH 


HV 


W HH HV W 


Boxcar 


15.696 


5.768 


25.111 0.083 0.038 0.083 


Refined Lee 


11.665 


10.136 


14.398 0.164 0.092 0.144 


IDAN 


2.164 


3.171 


1.977 0.199 0.137 0.188 


SDNLM 80% 


7.269 


5.999 


11.217 0.234 0.150 0.230 


SDNLM 90% 


8.786 


6.578 


13.559 0.181 0.101 0.177 


SDNLM 99% 


14.429 


7.129 


23.787 0.101 0.055 0.101 



The ENL was estimated on homogeneous areas far from edges, so no smudg- 
ing from other areas contaminated these values. As expected, the most intense 
blurring produces the best results with respect to this criterion: the Boxcar filter 
outperforms in two out of three bands and, when, it is not the best, the Refined 
Lee filter is. Regarding the Equivalent Number of Looks, IDAN produces worse 
results, but a good performance in SSIM index. The SDNLM filter improves 
with respect to the ENL criterion when the significance level increases. In order 
to make it competitive with the Boxcar, Refined Lee and IDAN filters, the ap- 
plication to real data was done using 77 = 90%. Regarding the SSIM index, our 
proposal consistently outperforms the other three filters and, as expected, the 
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smallest the significance of the test the better the performance is with respect 
to this criterion since the least the image is blurred. 

Physical-based simulation Sant'Anna et al. [23 proposed a methodology for 
simulating PolSAR imagery taking into account the electromagnetic character- 
istics of the targets and of the sensing system. The simulated images are more 
realistic than the ones obtained by merely stipulating the posterior distribu- 
tion given the classes, in particular spatial correlation among pixels emerges and 
mixture of classes in the borders as observed. 

Each simulated pixel is a complex scattering matrix based on a phantom 
image (an idealized cartoon model) with five distinct regions. Multifrequency 
sets of single-look PolSAR images have been generated in the L-, C- and X-bands, 
corresponding to 1.25, 5.3 and 9.6 GHz, respectively. The acquisition geometry 
is that of an airborne monostatic sensor flying at 6,000 m of altitude and 35 ° 
grazing angle imaging a 290 x 290 m^ area terrain. The 3.0 m spatial resolution 
and 2.8 m pixel spacing were set in the range and the azimuth directions. The 
data have 128 x 128 pixels; details in [27] . 

Figur e [5] sh ows the simulated images in the L-, C- and X-bands (Figures |5(a)[ 
|5(f)| and |5(k)[ resp.) and its filtered versions by the Boxcar (Figures |5(b)[ |5(g)] 
and [5(1)1 resp.). Refined Lee (Figures |5(c)[ |5(h)| and |5(m)[ resp ), ID AN (Fig- 
ures ]5td)| |5(i)] and |5(^ resp.) and SDNLM filters (Figures |5(^ |5(j)] and \5{o)\ 



resp.) with r] = 90%. 

The three main visual drawbacks of the Boxcar, Refined Lee and ID AN filters 
are noticeable in these results, namely the excessive blurring introduced by the 
former, and the pixellate effect produced by the latter. The SDNLM filter con- 
sistently presents a good compromise between smoothing and edge preservation. 

Table |2] presents the assessment of the filters in three intensity channels by 
means of the ENL (computed in the central region, a homogeneous area), and 
by means of the SSIM index. The best results are highlighted in boldface. The 
results are consistent with those observed in Table [ij the Boxcar filter is the best 
with respect to the ENL computed in homogeneous areas far from edges, but 
the SDNLM filter outperforms the other three filters when a more sophisticated 
metric is used: the SSIM index, which takes into account not only smoothness 
but structural information. Regarding these data, the ideal significance level lies 
approximately between r] = 80% and 77 = 90%; these values provide a good 
smoothing without compromising the structural information. 

6.2 Real data 

In the remainder of this section the data will be presented in false color using the 
Pauli decomposition [20]. This representation of PolSAR data has the advantage 
of being interpretable in terms of types of backscattering mechanisms. It consists 
of assigning \Shh-^Svv\'^ to the Red channel, \Shh — Svv\'^ to the Green channel, 
and 2\Shv\'^ to the Blue channel. This is one of many possible representations, 
and it is noteworthy that the data are filtered in their original domain, before 
being decomposed for visualization. 
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Table 2. Image quality indexes in the images shown in Figure [5] 



S Filter 




ENL 


SSIM Index 


Co 


HH 


HV 


W HH HV W 


Boxcar 


19.294 21.952 23.072 0.067 0.077 0.115 


Refined Lee 


12.639 


14.059 


14.850 0.159 0.163 0.186 


IDAN 


4.278 


4.400 


4.655 0.239 0.205 0.237 


^ SDNLM 80% 


8.683 


9.619 


9.588 0.246 0.231 0.251 


SDNLM 90% 


9.737 


11.128 


10.780 0.120 0.187 0.220 


SDNLM 99% 


17.559 


20.359 


21.125 0.100 0.168 0.140 


Boxcar 


24.008 21.970 23.547 0.078 0.077 0.087 


Refined Lee 


14.305 


12.947 


12.434 0.177 0.165 0.144 


IDAN 


2.626 


3.079 


4.492 0.268 0.253 0.243 


^ SDNLM 80% 


10.290 


9.217 


9.429 0.271 0.265 0.244 


SDNLM 90% 


11.840 


10.599 


10.801 0.212 0.211 0.194 


SDNLM 99% 


22.061 


19.576 


21.045 0.110 0.110 0.114 


Boxcar 


18.553 23.125 23.747 0.080 0.102 0.154 


Refined Lee 


9.694 


13.603 


14.526 0.169 0.195 0.235 


IDAN 
^ SDNLM 80% 


1.945 


2.151 


2.986 0.259 0.259 0.273 


8.348 


9.463 


9.574 0.267 0.270 0.300 


SDNLM 90% 


9.371 


11.049 


11.620 0.210 0.222 0.261 


SDNLM 99% 


17.023 


21.082 


21.547 0.121 0.131 0.190 



A National Aeronautics and Space Administration Jet Propulsion Laboratory 
(NASA/JPL) Airborne SAR (AIRSAR) image of the San Francisco Bay was used 
for evaluating the proposed filter, see |http : //earth . eo . esa . int/polsar pro^ 
'datasets-html, The original PolSAR image was generated in the L-band, four 
nominal looks, and 10 x 10 m spatial resolution. The test region has 350 x 350 
pixels, and is shown in Figure |6(b)| along with a Google Map© of the area 
(Figure [6(a)| see 'http : //gooTgl/maps/H JkPf ') . 

The four filters employ a kernel of 5 x 5 pixels, and the patches in our proposal 



are 3x3 windows. Figure 6(c) shows the effect of the Boxcar filter. Albeit the 
noise reduction is evident, it is also clear that the blurring introduced eliminates 
useful information as, for instance, the Presidio Golf Course: the curvilinear dark 
features to the center of the image: a forested area. Figure [6(d)] is the result of 
applying the Refined Lee filter, which shows a good performance, but some 
details in the edges are eliminated. In particular, the Mountain Lake, the small 
brown spot to the center of the image, is blurred, as well as the blocks in the 
urban area. The results of the IDAN filter and of our proposal with r] = 90% are 
shown in Figures |6(e)| and |6(f)[ respectively. Both filters are able to smooth the 
image in a selective way, but the SDNLM filter enhances more the signal-to-noise 
ratio while preserving fine details than the IDAN filter. 

In image quality assessment, the SSIM requires a reference image, as was the 



case of sections 6.1 and 6.1, therefore this index is not applied with ease on real 
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images. Table |3] presents the result of assessing the filters in the three intensity 
channels by means of the ENL in the large forest area, a homogeneous area, and 
of the BRISQUE index. Again, the results are consistent with what was observed 
before: the mere evaluation of the noise reduction by the ENL suggests the 
Boxcar filter as the best one, but the natural scene distortion-generic BRISQUE 
index is better after applying the SDNLM filter in all intensity channels. 

Table 3. Image quality indexes in the real PolSAR image. 



Filter 




ENL 


BRISQUE Index 




HH 


HV W 


HH HV W 


Real data 


3.867 


4.227 4.494 


58.258 70.845 61.593 


Boxcar 


14.564 25.611 18.946 36.498 37.714 36.792 


Refined Lee 


n.491 


20.415 15.407 


44.997 51.547 49.412 


IDAN 


2.994 


3.732 3.923 


28.823 34.853 34.691 


SDNLM 80% 


7.263 


11.532 8.299 


27.841 27.256 33.541 


SDNLM 90% 


8.177 


12.404 9.013 


28.622 35.622 35.016 


SDNLM 99% 


10.828 


18.379 13.075 


31.026 33.881 36.881 



6.3 Effect of the filters on the scattering characteristics 

Polarimetric target decompositions aim at expressing the physical properties 
of the scattering mechanisms. Among them, the entropy-based decomposition 
proposed by Cloude and Pottier [7 is widely used in PolSAR image classification. 
It extracts two parameters from each observed covariance matrix: the scattering 
entropy H G [0, 1], and a G [0,90], an indicator of the type of scattering. The 
{H^ a) plane is then divided into nine regions which provide both a classification 
rule and an interpretation for the observed data. 

Figure [T] shows the effect of the filters on the entropy of the data. Three classes 
of entropy are barely visible corresponding, in increasing brightness, to the sea, 
urban area an forest. The difference between the classes is small. All filters 
enhance the discrimination ability of the entropy, cf. figures Figures 7(b)[ |7(c)[ 



7(d) and|7(e)|with Figure 7(a) Again, the best preservation of detail is obtained 



by the IDAN and SDNLM techniques that, among others, retain the information 
of the Presidio Golf Course. Additionally, our proposal is the one that best 
preserves the low entropy small spots within the urban area, a feature typical of 
the high variability of these areas due to their heterogeneous composition. 

Samples from the sea, the urban area, and the forest, identified in red, blue 
and magenta, respectively, in Figure [6 (b)[ were taken. The (H^a) value of each 
point from these samples is presented in Figure [S] before and after applying the 
filters. The improvement of using filters is notorious when comparing Figure [8(a)] 
with figures [8(b)||8(c)[ 8(d)| and |8(e)[ While the original data is mixed, specially 
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the samples from urban and forest areas, after applying the filters the data tend 
to group in clusters. 

The sea samples are confined to Zone 9 in all datasets, and all filters have 
the effect of reducing their variability. While the Boxcar and Refined Lee filters 
produce very similar clusters of data, the ID AN and SDNLM filters reduce both 
the entropy and the a coefficient, but still within the zone of low entropy surface 
scatter, making the sample much more distinguishable from the rest of the data. 

The samples from urban (in blue) and forest (in magenta) areas have different 
mean values of entropy and a. The former occupy mostly zones 4 (medium 
entropy multiple scattering) and 5 (medium entropy vegetation scattering), while 
the latter span mostly zones 2 (high entropy vegetation scattering) and 5. While 
both are present in Zone 5, they seldom overlap; the forest samples have higher 
values of a. Comparing these two samples in the images filtered by the Refined 
Lee and SDNLM techniques, one notices that they would produce very similar 
classifications. The SDNLM produces clusters with more spread than the Refined 
Lee, but not at the expense of mixing different classes. 

In this manner, the filters preserve the scattering properties of the samples, a 
central feature of every speckle smoothing technique for PolS AR data, according 
to Lee and Pottier |20j. 

6.4 Effect of iterations number in filtering 

As previously discussed, SDNLM can be iterated since the properties upon which 
it is based are preserved by convolutions. Figure |9] presents the original image 
for reference (Figure [6 (b)[ ), and the result of applying each technique (Boxcar, 
Refined Lee, ID AN and SDNLM with 77 = 80% in each row) one, three and five 
times (first, second and third column, respectively). 

The most notorious new result stems from comparing the ID AN and SDNLM 
filters. Three iterations are enough for the former to smudge the original data, 
and with five iterations the blurring it produces is comparable with that of the 
Refined Lee and Boxcar filters. The SDNLM filter, even after five iterations, still 
preserves most of the spatial information. 



Figure 10 presents the {H, a) scatter plot of the samples before and after 
iterating the filters one, three and five times. Each new iteration adds cohesion to 
the clusters, whatever the filter employed. Regarding the SDNLM, the difference 
between one and three iterations is noticeable. 

Table [4] presents the quantitative analysis of the resulting images. Again, the 
Boxcar procedure yields the best noise reduction in smudge-free areas in most 
of the situations, followed by the Refined Lee filter. The SDNLM filter tuned 
somewhere between 80% and 90% of significance produces the best BRISQUE 
indexes. 

Table [5] presents the values of mean, variance and ENL estimator on the 
samples from the three regions of interest (sea, urban and forest), in the cross- 
polarized band (HV). The fist line in Table p] presents the values observed in 
the original (unfiltered) image. The best values are highlighted in bold; being 
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Table 4. Image quality indexes in the real PolSAR image. 



Filter 




ENL 




BRISQUE Index 




HH 


HV 


W 


HH HV W 


Real data 


3.867 


4.227 


4.494 


58.258 70.845 61.593 



Boxcar 
^0 Refined Lee 


14.564 25.611 18.946 


36.498 


37.714 


36.792 


11.491 


20.415 


15.407 


44.997 


51.547 


49.412 


1d IDAN 


2.994 


3.732 


3.923 


28.823 


34.853 


34.691 


1 SDNLM 80% 


7.263 


11.532 


8.299 


27.841 27.256 33.541 


'X SDNLM 90% 


8.177 


12.404 


9.013 


28.622 


35.622 


35.016 


SDNLM 99% 


10.828 


18.379 


13.075 


31.026 


33.881 


36.881 


y] Boxcar 


24.107 


54.348 


37.501 


78.379 


63.823 


79.307 


Refined Lee 


24.451 


45.088 


47.423 62.421 


42.882 


67.016 


•■g IDAN 


30.880 


40.220 


41.640 


55.811 


52.470 


67.016 


S SDNLM 80% 


15.502 


32.595 


21.686 


39.500 39.801 42.866 


•T SDNLM 90% 


17.279 


35.896 


23.478 


41.800 


42.371 


45.177 


^ SDNLM 99% 


21.892 


50.069 


33.347 


55.824 


52.824 


56.844 


y] Boxcar 


29.950 


73.504 50.811 


81.585 


71.628 


82.365 


Refined Lee 


36.142 


72.948 


91.634 70.908 


51.728 


77.271 


^ IDAN 


31.310 


41.050 


42.340 


65.906 


64.803 


66.400 


S SDNLM 80% 


20.069 


52.620 


31.233 


49.357 


46.318 


50.872 


•T SDNLM 90% 


23.793 


61.817 


36.548 


48.803 47.177 


50.430 


"^ SDNLM 99% 


27.798 


72.209 


46.450 


66.251 


58.663 


65.944 



the best mean the closest to the original value, and the best standard deviation 
the smallest one. 

The SDNLM filter is the best at preserving the original mean values, and 
with reduced standard deviation. Our proposal does not provide the best vari- 
ability reduction, a behavior that may be associated with the preservation of fine 
structures and a smaller loss of spatial resolution, as can be noted in Figure [9J 
The equivalent number of looks behaves consistently with what was observed in 
previous examples. 

7 Conclusions 

The use h-(j) divergences, a tool in Information Theory, led to test statistics (with 
a known and tractable asymptotic distribution) able to check if two samples 
cannot be described by the same complex Wishart distribution, the classical 
model for PolSAR data. Using one of these test statistics, namely the one based 
on the Hellinger distance, we devised a convolution filter whose weights are 
function of the p-value of tests which compare two patches of size 3 x 3 in a 
search window of size 5x5 pixels. 

The filter obtained in this manner (SDNLM - Stochastic Distances Nonlocal 
Means) was compared with the Boxcar, Refined Lee and IDAN filters in a variety 
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of PolS AR imagery: data simulated from the complex Wishart law over a realistic 
phantom using parameters observed in practice, data simulated from the elec- 
tromagnetic properties of the scattering over a simplified cartoon model, and a 
real PolSAR image over San Francisco, CA. The quantitative assessment verified 
the equivalent number of looks (a measure of noise reduction) over smudge-free 
samples, the structural SSIM index, and the BRISQUE index used appropriately 
on no-references images (real images or blind). The Boxcar filter promotes the 
strongest noise reduction in these conditions, but at the expense of obliterating 
small details. The Refined Lee and IDAN filters are competitive, but produce a 
pixellated effect and their SSIM index is worse than the produced by the SDNLM 
filter in all instances. We noted this same feature with the BRISQUE index ap- 
plied to real data and the index values remain stable even during iteration of 
the SDNLM filter, which does not happen with other filters assessed. 

A qualitative assessment was also made checking how the polarimetric en- 
tropy is affected by the filters. We noticed that all the filters enhance it but, 
in particular, our proposal performs the most refined enhancement since it pre- 
serves very small details which are characteristic of complex urban areas. 

The effect of the filters and of applying them iteratively was also verified 
in the (H^a) plane. All filters produced more and more compact clusters of 
observations in this plane as more iterations were applied. The SDNLM filter 
yielded the best separation of the sea sample, while the other two were treated 
at least as well as they were by the other filters. 

The SDNLM filter has three tuning parameter: (i) filtering window size, 
(ii) the size of the patches, and (iii) the significance of the test. We provide a 
range of suggested values for the latter, and show good results with an economic 
choice for the two former. The filter can also be applied iteratively if more 
smoothing is requested and, provided an adequate statistical model, it can also 
be applied to other types of data. 

We conclude that our proposal is a good candidate for smoothing PolSAR im- 
agery without compromising either small details or the scattering characteristics 
of the targets. 

The test statistics are invariant with respect to permutations of the sample; 
directional features will be considered in forthcoming works. Future research 
includes the proposal of quality measures for PolSAR imagery, and the use of 
tests based on entropies [12j. 
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A Computational information 



Computing the ML estimator of the equivalent number of looks L given in equa- 
tion ([2| and the test statistic based on the Hellinger distance, c.f. equation (p|, 
are part of the computational core of this proposal. Each weight requires com- 
puting these quantities. It is noteworthy that they involve only two operations 
on complex matrices: the determinant and the inverse. Using the fact that the 
matrices are Hermitian and positive definite, it is possible to reduce drastically 
the number of operations required to calculate these two quantities. Specialized 
accelerated function in the R programming language [31 were developed with 
the speed of the filters in mind. 

The time required to filter a 128 x 128 pixels image with one iteration is 
of about 75 s in an Intel® Core^^ i7-3632QM CPU 2.20 GHz, with software 
developed in R version 2.14.1 running on Ubuntu 12.04. R was the choice because 
of its excellent accuracy with respect to similar platforms pLj . 



B Observed covariance matrices 



In the following we present the observed covariance matrices that were used to 
simulate the data presented in Figure |3(b)[ 
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^1 



7.60830 -0.74901 - 2.29165i 1.38157 + 8.39200i 
24.8580 -5.90346 - 0.45011i 
32.2771 

128.592 12.1941 - 7.12246i 39.1107 + 18.7954i " 
336.959 -8.49716 - 11.8210i 
154.343 

29.6303 4.86985 + 1.55848i 3.41851 + 1.43502i " 
86.8985 -2.03628 - 8.24319i 
43.3504 

14.0576 -0.25731 - 1.48967i 4.36926 + 9.41493i 
60.5614 -4.92951 - 2.16850i 
42.3767 

4.89301 -0.52225 - 0.62765i 1.38866 + 5.29889i 
12.1149 -3.30897 - 0.85846i 
25.6761 

18.7013 0.81235 - 1.72513Z 1.26677 + 6.08878i 
32.8094 -3.01618 - 1.67916i 
25.8651 
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(a) Phantom (b) Simulated single-look data 

^1 





(c) Boxcar filter 



(d) Refined Lee filter 




(e) IDAN filter 



(f) SDNLM filter 



Fig. 3. Original single-look PolSAR simulated data and filtered versions after one it- 
eration. 
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(a) Phantom edges 



(b) Edges of the simulated HH band 




(c) Boxcar filter edges 



(d) Refined Lee filter edges 




(e) IDAN filter edges 



(f ) SDNLM filter edges 



Fig. 4. Edges detected by the Canny filter applied to the HH polarization channel of 
the original and filtered images. 
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(a) Simulated (b) Boxcar fil- (c) Refined Lee (d) IDAN filter (e) SDNLM fil- 
L-band ter filter ter 



l^'Tfi^'^f^ 




(f) Simulated (g) Boxcar fil- (h) Refined Lee (i) IDAN filter (j) SDNLM fil- 
C-band ter filter ter 




(k) Simulated (1) Boxcar filter (m) Refined (n) IDAN filter (o) SDNLM fil- 
X-band Lee filter ter 



Fig. 5. Original physical-based single- look images and their filtered versions after one 
iteration. 
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(a) Map of the area 



(b) AIRSAR L-band data 




(e) IDAN filter 



(f) SDNLM filter 



Fig. 6. Pauli decomposition of the original AIRSAR image over San Francisco and its 
filtered versions. 
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(a) Entropy of the Original data 




(b) Entropy of the Boxcar filter (c) Entropy of the Refined Lee filter 




(d) Entropy of the IDAN filter (e) Entropy of the SDNLM filter 

Fig. 7. Entropy in the AIRSAR L-band data. 
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(b) Boxcar filtered data 



(c) Refined Lee filtered data 
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(d) ID AN filtered data 
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(e) SDNLM filtered data 



Fig. 8. Scatter plot in the (i/, a) plane of samples from the AIRSAR L-band image. 
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(a) Boxcar filter 7=1 (b) Boxcar filter 7 = 3 (c) Boxcar filter 7 = 5 




(d) Refined Lee filter 7 = 1 (e) Refined Lee filter 7 = 3 (f ) Refined Lee filter 1 — 5 




(g) IDAN filter 1=1 (h) IDAN filter 1=1 (i) IDAN filter 7 = 1 




(j) SDNLM filter 1=1 (k) SDNLM filter 7 = 3 (1) SDNLM filter 7 = 5 



Fig. 9. Pauli decomposition of images filtered one, three and five times with each 
technique. 
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(a) AIRSAR L-band data 
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(b) Boxcar filter 7=1 (c) Boxcar filter 7 = 3 (d) Boxcar filter 7 = 5 
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(e) Refined Lee filter 7 = 1 (f ) Refined Lee filter 7 = 3 (g) Refined Lee filter 1 — 5 
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(h) ID AN filter 7=1 (i) ID AN filter 7 = 3 (j) ID AN filter 7 = 5 
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(k) SDNLM filter 7=1 (1) SDNLM filter 7 = 3 (m) SDNLM filter 7 = 5 



Fig. 10. Scatter plot in the (77, a) plane of samples from the AIRSAR L-band image 
after one, three and five iterations 7 
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